knitr::opts_chunk$set(echo = TRUE)

#cmdstanr::set_cmdstan_path(path = "C:/Users/kueng/.cmdstan/cmdstan-2.35.0")
#cmdstanr::set_cmdstan_path(path = "C:/Users/pascku/.cmdstan/cmdstan-2.36.0")

library(tidyverse)
library(R.utils)
library(wbCorr)
library(readxl)
library(kableExtra)
library(brms)
library(bayesplot)
library(see)
library(beepr)
library(DHARMa)
library(digest)



source(file.path('Functions', 'ReportModels.R'))
source(file.path('Functions', 'PrettyTables.R'))
source(file.path('Functions', 'ReportMeasures.R'))
source(file.path('Functions', 'PrepareData.R'))

report_function_hash <- digest::digest(summarize_brms)
system("shutdown /a")
## [1] 1116
# Set options for analysis
use_mi = FALSE
shutdown = FALSE
report_ordinal = FALSE
do_priorsense = FALSE
get_bayesfactor = TRUE
check_models = TRUE # ! set to TRUE !

if (get_bayesfactor) {
  stats_to_report <- c('CI', 'SE', 'pd', 'ROPE', 'BF', 'Rhat', 'ESS')
} else {
  stats_to_report <- c('CI', 'SE', 'pd', 'ROPE', 'Rhat', 'ESS')
}

options(
  dplyr.print_max = 100, 
  brms.backend = 'cmdstan',
  brms.file_refit = ifelse(use_mi, 'never', 'on_change'),
  brms.file_refit = 'on_change',
  #brms.file_refit = 'always',
  error = function() {
    beepr::beep(sound = 5)
    if (shutdown) {
      system("shutdown /s /t 180")
      quit(save = "no", status = 1)
    }
  }
  , es.use_symbols = TRUE
)


####################### Model parameters #######################

iterations = 12000 # 12'000 per chain to achieve 40'000
warmup = 2000 # 2000

# NO AR!!!
#corstr = 'ar'
#corstr = 'cosy_couple'
#corstr = 'cosy_couple:user'


################################################################

suffix = paste0('_final_with_plan_', as.character(iterations))
df <- openxlsx::read.xlsx(file.path('long.xlsx'))
df_original <- df

df_double <- prepare_data(df, recode_pushing = TRUE, use_mi = use_mi)[[1]]

Constructing scales Re-coding pusing reshaping data (4field) centering data within and between

summary(df_double$pushing)

Min. 1st Qu. Median Mean 3rd Qu. Max. NA’s 0.0000 0.0000 0.0000 0.1635 0.0000 5.0000 241

Modelling

# For indistinguishable Dyads
model_rows_fixed <- c(
    'Intercept', 
    # '-- WITHIN PERSON MAIN EFFECTS --', 
    'persuasion_self_cw', 
    'persuasion_partner_cw', 
    'pressure_self_cw', 
    'pressure_partner_cw', 
    'pushing_self_cw', 
    'pushing_partner_cw', 
    'day', 
    'plan_selfPlan',
    'plan_partnerPlan',
    'weartime_self_cw',
    
    # '-- BETWEEN PERSON MAIN EFFECTS',
    'persuasion_self_cb',
    'persuasion_partner_cb',
    'pressure_self_cb',
    'pressure_partner_cb',
    'pushing_self_cb',
    'pushing_partner_cb',
    'weartime_self_cb'
  )


model_rows_fixed_ordinal <- c(
  model_rows_fixed[1],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rows_fixed[2:length(model_rows_fixed)]
)

model_rows_random <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  'sd(persuasion_self_cw)',
  'sd(persuasion_partner_cw)',
  'sd(pressure_self_cw)',
  'sd(pressure_partner_cw)',
  'sd(pushing_self_cw)',
  'sd(pushing_partner_cw)',
  # '-- CORRELATION STRUCTURE -- ', 
  'sigma'
)

model_rows_random_ordinal <- c(model_rows_random,'disc')
# For indistinguishable Dyads
model_rownames_fixed <- c(
    "Intercept", 
    # "-- WITHIN PERSON MAIN EFFECTS --", 
    "Daily persuasion experienced", 
    "Daily persuasion utilized (partner's view)", # OR partner received
    "Daily pressure experienced", 
    "Daily pressure utilized (partner's view)", 
    "Daily pushing experienced", 
    "Daily pushing utilized (partner's view)", 
    "Day", 
    "Own Actionplan",
    'Partner Actionplan',
    "Daily weartime",
    
    # "-- BETWEEN PERSON MAIN EFFECTS",
    "Mean persuasion experienced", 
    "Mean persuasion utilized (partner's view)", 
    "Mean pressure experienced", 
    "Mean pressure utilized (partner's view)", 
    "Mean pushing experienced", 
    "Mean pushing utilized (partner's view)", 
    "Mean weartime"
  )


model_rownames_fixed_ordinal <- c(
  model_rownames_fixed[1],
  'Intercept[1]',
  'Intercept[2]',
  'Intercept[3]',
  'Intercept[4]',
  'Intercept[5]',
  model_rownames_fixed[2:length(model_rownames_fixed)]
)

model_rownames_random <- c(
  # '--------------',
  # '-- RANDOM EFFECTS --',
  'sd(Intercept)', 
  "sd(Daily persuasion experienced)", 
  "sd(Daily persuasion utilized (partner's view))", # OR partner received
  "sd(Daily pressure experienced)", 
  "sd(Daily pressure utilized (partner's view))", 
  "sd(Daily pushing experienced)", 
  "sd(Daily pushing utilized (partner's view))", 
  # '-- CORRELATION STRUCTURE -- ', 
  'sigma'
)

model_rownames_random_ordinal <- c(model_rownames_random,'disc')
rows_to_pack <- list(
  "Within-Person Effects" = c(2,11),
  "Between-Person Effects" = c(12,18),
  "Random Effects" = c(19, 25), 
  "Additional Parameters" = c(26,26)
  )


rows_to_pack_ordinal <- list(
  "Within-Person Effects" = c(2+5,11+5),
  "Between-Person Effects" = c(12+5,18+5),
  "Random Effects" = c(19+5, 25+5), 
  "Additional Parameters" = c(26+5,26+6)
  )

Self-Reported MVPA

range(df_double$pa_sub, na.rm = T) 
## [1]   0 720
hist(df_double$pa_sub, breaks = 40) 

hist(log(df_double$pa_sub+00000000001), breaks = 40)

Hurdle Lognormal Model

formula <- bf(
  pa_sub ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    plan_self + plan_partner +
    day + 
    
    # Random effects
    (1 + persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | dd | coupleID),
  
  hu = ~ persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    plan_self + plan_partner +
    day + 
    
    # Random effects
    (1 + persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | dd | coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
) 

prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 2)", class = "b", dpar = "hu")
  , brms::set_prior("normal(0, 50)", class = "Intercept") # for non-zero PA
  , brms::set_prior("normal(0.5, 2.5)", class = "Intercept", dpar = 'hu') # hurdle part
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = hurdle_lognormal()
#)

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_sub <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::hurdle_lognormal(), 
  #family = brms::hurdle_negbinomial(), 
  #family = brms::hurdle_poisson(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 42,
  file = file.path("models_cache_brms", paste0("pa_sub_hu_lognormal", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
pa_sub_digest <- digest::digest(pa_sub)
if (check_models) {
  check_brms(pa_sub, log_pp_check = TRUE)
  DHARMa.check_brms.all(pa_sub, integer = TRUE, outliers_type = 'bootstrap')
}
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                   Term  VIF   VIF 95% CI Increased SE Tolerance
##     persuasion_self_cw 1.56 [1.51, 1.62]         1.25      0.64
##  persuasion_partner_cw 1.08 [1.06, 1.12]         1.04      0.92
##       pressure_self_cw 1.06 [1.04, 1.10]         1.03      0.94
##    pressure_partner_cw 1.06 [1.04, 1.10]         1.03      0.94
##        pushing_self_cw 1.06 [1.04, 1.10]         1.03      0.94
##     pushing_partner_cw 1.08 [1.05, 1.12]         1.04      0.93
##     persuasion_self_cb 1.09 [1.06, 1.12]         1.04      0.92
##  persuasion_partner_cb 4.70 [4.48, 4.93]         2.17      0.21
##       pressure_self_cb 4.60 [4.39, 4.83]         2.15      0.22
##    pressure_partner_cb 3.26 [3.12, 3.41]         1.81      0.31
##        pushing_self_cb 3.16 [3.02, 3.30]         1.78      0.32
##     pushing_partner_cb 2.81 [2.69, 2.93]         1.68      0.36
##              plan_self 2.73 [2.61, 2.85]         1.65      0.37
##           plan_partner 1.37 [1.33, 1.42]         1.17      0.73
##                    day 1.37 [1.32, 1.42]         1.17      0.73
##  Tolerance 95% CI
##      [0.62, 0.66]
##      [0.89, 0.94]
##      [0.91, 0.96]
##      [0.91, 0.96]
##      [0.91, 0.96]
##      [0.90, 0.95]
##      [0.89, 0.94]
##      [0.20, 0.22]
##      [0.21, 0.23]
##      [0.29, 0.32]
##      [0.30, 0.33]
##      [0.34, 0.37]
##      [0.35, 0.38]
##      [0.70, 0.75]
##      [0.71, 0.76]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy        100%
## 
## Predicted Distribution of Response
## 
##                Distribution Probability
##  neg. binomial (zero-infl.)         78%
##               beta-binomial         16%
##                   lognormal          3%
## 
## Divergences:
## 0 of 40000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 40000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 1 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 7, observations = 3736, p-value = 0.8
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.000267666 0.002944325
## sample estimates:
## outlier frequency (expected: 0.00157922912205567 ) 
##                                        0.001873662
if (do_priorsense) {
  priorsense_vars <- c(
      'Intercept',
      'b_persuasion_self_cw',
      'b_persuasion_partner_cw',
      'b_pressure_self_cw',
      'b_pressure_partner_cw',
      'b_pushing_self_cw',
      'b_pushing_partner_cw'
  )
  
  hurdle_priorsense_vars <- c(
    'Intercept_hu',
    'b_hu_persuasion_self_cw',
    'b_hu_persuasion_partner_cw',
    'b_hu_pressure_self_cw',
    'b_hu_pressure_partner_cw',
    'b_hu_pushing_self_cw',
    'b_hu_pushing_partner_cw'
  )
  
  gc()
  priorsense::powerscale_sensitivity(pa_sub, variable = c(priorsense_vars, hurdle_priorsense_vars))
  priorsense::powerscale_plot_dens(pa_sub, variable = c(priorsense_vars, hurdle_priorsense_vars))
  priorsense::powerscale_plot_ecdf(pa_sub, variable = c(priorsense_vars, hurdle_priorsense_vars))
  priorsense::powerscale_plot_quantities(pa_sub, variable = c(priorsense_vars, hurdle_priorsense_vars))
}
# rope range for continuous part of the model
rope_factor <- sd(log(pa_sub$data$pa_sub[pa_sub$data$pa_sub > 0]))
rope_range_continuous = c(-0.1 * rope_factor, 0.1 * rope_factor)

summary_pa_sub <- summarize_brms(
  pa_sub, 
  stats_to_report = stats_to_report,
  rope_range = rope_range_continuous,
  hu_rope_range = c(-0.18, 0.18),
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = T) 
## Sampling priors, please wait...
## Warning in summarize_brms(pa_sub, stats_to_report = stats_to_report, rope_range
## = rope_range_continuous, : Coefficients were exponentiated. Double check if
## this was intended.
# Print the updated dataframe
summary_pa_sub %>%
  print_df(rows_to_pack = rows_to_pack)
exp(Est.)_hu SE_hu 95% CI_hu pd_hu ROPE_hu inside ROPE_hu BF_hu BF_Evidence_hu Rhat_hu Bulk_ESS_hu Tail_ESS_hu exp(Est.)_nonzero SE_nonzero 95% CI_nonzero pd_nonzero ROPE_nonzero inside ROPE_nonzero BF_nonzero BF_Evidence_nonzero Rhat_nonzero Bulk_ESS_nonzero Tail_ESS_nonzero
Intercept 0.26*** 0.04 [ 0.19, 0.35] 1.000 [0.84, 1.20] 0.000 >100 Overwhelming Evidence 1.000 17636 25222 37.84*** 2.60 [33.00, 43.33] 1.000 [0.92, 1.08] 0.000 >100 Overwhelming Evidence 1.001 9733 17297
Within-Person Effects
Daily persuasion experienced 1.58*** 0.11 [ 1.39, 1.83] 1.000 [0.84, 1.20] 0.000 >100 Overwhelming Evidence 1.000 25310 25796 1.03 0.03 [ 0.98, 1.09] 0.868 [0.92, 1.08] 0.963 0.017 Very Strong Evidence for Null 1.000 13834 23907
Daily persuasion utilized (partner’s view) 1.34*** 0.08 [ 1.19, 1.52] 1.000 [0.84, 1.20] 0.035 >100 Overwhelming Evidence 1.000 27878 24532 1.03 0.02 [ 0.99, 1.08] 0.928 [0.92, 1.08] 0.974 0.022 Very Strong Evidence for Null 1.000 19254 25751
Daily pressure experienced 0.96 0.15 [ 0.69, 1.30] 0.614 [0.84, 1.20] 0.739 0.079 Strong Evidence for Null 1.000 27614 23572 0.91 0.05 [ 0.82, 1.00] 0.972 [0.92, 1.08] 0.345 0.051 Strong Evidence for Null 1.000 31401 26663
Daily pressure utilized (partner’s view) 1.49* 0.28 [ 1.05, 2.34] 0.986 [0.84, 1.20] 0.111 1.196 Weak Evidence 1.000 26732 19474 0.95 0.04 [ 0.86, 1.03] 0.892 [0.92, 1.08] 0.700 0.014 Very Strong Evidence for Null 1.000 30682 28039
Daily pushing experienced 0.95 0.14 [ 0.71, 1.30] 0.636 [0.84, 1.20] 0.740 0.079 Strong Evidence for Null 1.000 20946 22281 0.99 0.03 [ 0.92, 1.05] 0.640 [0.92, 1.08] 0.969 0.008 Very Strong Evidence for Null 1.000 26058 25658
Daily pushing utilized (partner’s view) 1.29* 0.14 [ 1.05, 1.62] 0.993 [0.84, 1.20] 0.247 1.023 Weak Evidence 1.000 29406 25834 0.96 0.03 [ 0.90, 1.02] 0.916 [0.92, 1.08] 0.878 0.018 Very Strong Evidence for Null 1.000 27923 27480
Day 0.86 0.13 [ 0.64, 1.16] 0.834 [0.84, 1.20] 0.576 0.118 Moderate Evidence for Null 1.000 45401 29525 0.99 0.06 [ 0.88, 1.12] 0.577 [0.92, 1.08] 0.786 0.007 Very Strong Evidence for Null 1.000 39959 29261
Own Actionplan 9.45*** 0.96 [ 7.76, 11.57] 1.000 [0.84, 1.20] 0.000 >100 Overwhelming Evidence 1.000 40464 31867 1.32*** 0.06 [ 1.21, 1.45] 1.000 [0.92, 1.08] 0.000 >100 Overwhelming Evidence 1.000 43927 29958
Partner Actionplan 1.17 0.12 [ 0.96, 1.42] 0.940 [0.84, 1.20] 0.601 0.167 Moderate Evidence for Null 1.000 37878 31193 1.08 0.05 [ 0.99, 1.17] 0.956 [0.92, 1.08] 0.541 0.031 Strong Evidence for Null 1.000 37238 29764
Daily weartime NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 1.49 0.54 [ 0.72, 3.06] 0.863 [0.84, 1.20] 0.216 0.332 Weak Evidence for Null 1.000 9680 17445 0.99 0.15 [ 0.74, 1.34] 0.514 [0.92, 1.08] 0.394 0.014 Very Strong Evidence for Null 1.000 7913 15602
Mean persuasion utilized (partner’s view) 1.40 0.50 [ 0.69, 2.85] 0.826 [0.84, 1.20] 0.252 0.280 Moderate Evidence for Null 1.000 9860 17614 0.95 0.14 [ 0.71, 1.28] 0.636 [0.92, 1.08] 0.371 0.015 Very Strong Evidence for Null 1.000 7546 15092
Mean pressure experienced 0.39* 0.16 [ 0.17, 0.88] 0.988 [0.84, 1.20] 0.030 2.695 Weak Evidence 1.000 14859 22272 1.21 0.21 [ 0.86, 1.69] 0.861 [0.92, 1.08] 0.203 0.019 Very Strong Evidence for Null 1.000 12502 22414
Mean pressure utilized (partner’s view) 0.49 0.21 [ 0.21, 1.14] 0.951 [0.84, 1.20] 0.087 0.845 Weak Evidence for Null 1.000 15616 22856 0.94 0.17 [ 0.66, 1.33] 0.638 [0.92, 1.08] 0.318 0.012 Very Strong Evidence for Null 1.000 10957 20414
Mean pushing experienced 1.10 0.59 [ 0.38, 3.17] 0.573 [0.84, 1.20] 0.261 0.269 Moderate Evidence for Null 1.000 13848 21535 1.24 0.28 [ 0.80, 1.94] 0.838 [0.92, 1.08] 0.174 0.019 Very Strong Evidence for Null 1.001 9386 17470
Mean pushing utilized (partner’s view) 1.96 1.04 [ 0.69, 5.73] 0.899 [0.84, 1.20] 0.120 0.586 Weak Evidence for Null 1.000 13955 21360 1.34 0.31 [ 0.85, 2.11] 0.903 [0.92, 1.08] 0.118 0.031 Strong Evidence for Null 1.000 8966 16483
Mean weartime NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.74 0.10 [0.57, 0.99] NA NA NA NA NA 1.000 16046 21340 0.30 0.04 [0.23, 0.40] NA NA NA NA NA 1.000 9362 16350
sd(Daily persuasion experienced) 0.21 0.09 [0.03, 0.39] NA NA NA NA NA 1.001 10442 7621 0.11 0.02 [0.07, 0.17] NA NA NA NA NA 1.000 17487 25737
sd(Daily persuasion utilized (partner’s view)) 0.15 0.09 [0.01, 0.34] NA NA NA NA NA 1.000 9486 10913 0.08 0.02 [0.04, 0.13] NA NA NA NA NA 1.000 11687 7419
sd(Daily pressure experienced) 0.20 0.18 [0.01, 0.74] NA NA NA NA NA 1.000 15688 17682 0.07 0.06 [0.00, 0.23] NA NA NA NA NA 1.000 13071 17615
sd(Daily pressure utilized (partner’s view)) 0.24 0.22 [0.01, 0.94] NA NA NA NA NA 1.000 14003 16516 0.06 0.05 [0.00, 0.18] NA NA NA NA NA 1.000 16359 14468
sd(Daily pushing experienced) 0.59 0.17 [0.31, 1.01] NA NA NA NA NA 1.000 16943 23757 0.09 0.04 [0.01, 0.17] NA NA NA NA NA 1.000 9282 7265
sd(Daily pushing utilized (partner’s view)) 0.18 0.14 [0.01, 0.50] NA NA NA NA NA 1.000 14479 13486 0.07 0.04 [0.01, 0.15] NA NA NA NA NA 1.000 11065 10435
Additional Parameters
sigma NA NA NA NA NA NA NA NA NA NA NA 0.68 0.01 [0.65, 0.70] NA NA NA NA NA 1.000 42655 29595 
# Plot continuous part of model

variable <- c(
  '(Intercept)',
  'b_persuasion_self_cw',
  'b_persuasion_partner_cw',
  'b_pressure_self_cw',
  'b_pressure_partner_cw',
  'b_pushing_self_cw',
  'b_pushing_partner_cw'
)


plot(
  bayestestR::p_direction(pa_sub, parameter = variable),
  priors = TRUE
) + theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(
    pa_sub, 
    parameter = variable, 
    range = rope_range_continuous,
    verbose = F,
    ci = 1
  )
) + theme_bw()

# Hurdle part of the model
variable <- c(
  'b_hu_persuasion_self_cw',
  'b_hu_persuasion_partner_cw',
  'b_hu_pressure_self_cw',
  'b_hu_pressure_partner_cw',
  'b_hu_pushing_self_cw',
  'b_hu_pushing_partner_cw'
)

plot(
  bayestestR::p_direction(pa_sub, parameter = variable),
  priors = TRUE
) + theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length

# The rope range for the bernoulli part of the model is -0.18, 0.18
plot(
  bayestestR::rope(pa_sub, parameter = variable, range = c(-0.18, 0.18), ci = 1),
  verbose = FALSE
) + theme_bw()
## Possible multicollinearity between b_hu_persuasion_partner_cb and
##   b_hu_persuasion_self_cb (r = 0.75), b_persuasion_partner_cb and
##   b_persuasion_self_cb (r = 0.75). This might lead to inappropriate
##   results. See 'Details' in '?rope'.

Hurdle part of the model on the left, non-zero part towards the right side of the table

conds_eff <- conditional_spaghetti(
  pa_sub, 
  effects = c(
    'persuasion_self_cw',
    'persuasion_partner_cw',
    'pressure_self_cw',
    'pressure_partner_cw',
    'pushing_self_cw',
    'pushing_partner_cw'
  ),
  x_label = c(
    'Received Persuasion',
    'Exerted Persuasion',
    'Received Pressure',
    'Exerted Pressure',
    'Received Plan-Related Pushing',
    'Exerted Plan-Related Pushing'
  ),
  group_var = 'coupleID',
  plot_full_range = TRUE,
  y_limits = c(0, 100),
  y_label = "Same-Day MVPA",
  y_labels = c('Probability of Being Active', 'Minutes of MVPA When Active', 'Overall Expected Minutes of MVPA'),
  , filter_quantiles = .9995
  , font_family = 'Candara'
)
## This is posterior version 1.6.0
## 
## Attaching package: 'posterior'
## The following object is masked from 'package:bayesplot':
## 
##     rhat
## The following objects are masked from 'package:stats':
## 
##     mad, sd, var
## The following objects are masked from 'package:base':
## 
##     %in%, match
## Registering fonts with R
## Warning: Dropping 'draws_df' class as required metadata was removed.
## Warning: Dropping 'draws_df' class as required metadata was removed.
## Warning: Dropping 'draws_df' class as required metadata was removed.
## Warning: Dropping 'draws_df' class as required metadata was removed.
## Warning: Dropping 'draws_df' class as required metadata was removed.
## Warning: Dropping 'draws_df' class as required metadata was removed.
print(conds_eff)

$persuasion_self_cw

## Warning in grid.Call(C_stringMetric, as.graphicsAnnot(x$label)): font family
## not found in Windows font database
## Warning in grid.Call(C_stringMetric, as.graphicsAnnot(x$label)): font family
## not found in Windows font database
## Warning in grid.Call(C_stringMetric, as.graphicsAnnot(x$label)): font family
## not found in Windows font database
## Warning: Removed 73 rows containing missing values or values outside the scale range
## (`geom_line()`).
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning: Removed 27 rows containing missing values or values outside the scale range
## (`geom_line()`).
## Picking joint bandwidth of 0.00809
## Warning in grid.Call(C_stringMetric, as.graphicsAnnot(x$label)): font family
## not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_stringMetric, as.graphicsAnnot(x$label)): font family
## not found in Windows font database
## Warning in grid.Call(C_stringMetric, as.graphicsAnnot(x$label)): font family
## not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database

$persuasion_partner_cw

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning: Removed 51 rows containing missing values or values outside the scale range
## (`geom_line()`).
## Warning: Removed 7 rows containing missing values or values outside the scale range
## (`geom_line()`).
## Picking joint bandwidth of 0.0063
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database

$pressure_self_cw

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning: Removed 1 row containing missing values or values outside the scale range
## (`geom_line()`).
## Picking joint bandwidth of 0.0112
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database

$pressure_partner_cw

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Picking joint bandwidth of 0.0185
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database

$pushing_self_cw

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Picking joint bandwidth of 0.0106
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database

$pushing_partner_cw

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Picking joint bandwidth of 0.00979
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database
## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, :
## font family not found in Windows font database

Note. This graphic illustrates the relationship between social control and moderate to vigorous physical activity (MVPA) using a Bayesian Hurdle-Lognormal Multilevel Model. The predictor is centered within individuals to examine how deviations from their average social control relate to same-day MVPA. Shaded areas indicate credible intervals, thick lines show fixed effects, and thin lines represent random effects, highlighting variability across couples. The plots display the probability of being active, expected minutes of MVPA when active, and combined predicted MVPA. The bottom density plot visualizes the posterior distributions of slope estimates, transformed to represent multiplicative changes in odds ratios (hurdle component) or expected values. Medians and 95% credible intervals (2.5th and 97.5th percentiles) are shown. Effects are significant, when the 95% credible interval does not overlap 1.

x_label = c(
    'Received Persuasion',
    'Exerted Persuasion',
    'Received Pressure',
    'Exerted Pressure',
    'Received Plan-Related Pushing',
    'Exerted Plan-Related Pushing'
  )

home_dir <- getwd()
output_dir <- file.path(home_dir, 'Output', 'Plots')

for (i in 1:length(conds_eff)) {
  effname <- names(conds_eff)[i]
  eff_plot <- conds_eff[[i]]
  x_label_i <- x_label[[i]]
  
  rmarkdown::render(
    file.path(output_dir, 'BeautifulPlotWithNote.Rmd'), 
    output_file = file.path(output_dir, paste0('Graphic_', effname, '.pdf')),
    params = list(
      home_dir = home_dir,
      output_dir = output_dir,
      p_i = eff_plot,
      p_name = effname,
      x_label = x_label_i
      ),
    envir = new.env(),
    quiet = TRUE
  )
}

print('done')

Comparing effect size of pressure and pushing

hypothesis(pa_sub, "pressure_self_cw < pushing_self_cw")
## Hypothesis Tests for class b:
##                 Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (pressure_self_cw... < 0    -0.09      0.07    -0.19     0.02        9.8
##   Post.Prob Star
## 1      0.91     
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.

Device Based MVPA

range(df_double$pa_obj, na.rm = T) 
## [1]   5.75 971.25
hist(df_double$pa_obj, breaks = 50)

df_double$pa_obj_log <- log(df_double$pa_obj)

hist(df_double$pa_obj_log, breaks = 50)

Lognormal Model

formula <- bf(
  pa_obj ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    plan_self + plan_partner +
    day + weartime_self_cw + weartime_self_cb +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
)



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 50)", class = "Intercept") 
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = lognormal()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

pa_obj_log <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = lognormal(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("pa_obj_log_gaussian", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
pa_obj_log_digest <- digest::digest(pa_obj_log)
if (check_models) {
  check_brms(pa_obj_log, log_pp_check = TRUE)
  DHARMa.check_brms.all(pa_obj_log, integer = TRUE, outliers_type = 'bootstrap')
}
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                   Term  VIF    VIF 95% CI Increased SE Tolerance
##     persuasion_self_cw 1.11 [1.08,  1.16]         1.05      0.90
##  persuasion_partner_cw 1.12 [1.09,  1.17]         1.06      0.89
##       pressure_self_cw 1.08 [1.05,  1.12]         1.04      0.93
##    pressure_partner_cw 1.11 [1.07,  1.15]         1.05      0.90
##        pushing_self_cw 1.13 [1.10,  1.18]         1.06      0.88
##     pushing_partner_cw 1.18 [1.14,  1.23]         1.09      0.85
##       pressure_self_cb 4.03 [3.82,  4.25]         2.01      0.25
##    pressure_partner_cb 4.34 [4.12,  4.59]         2.08      0.23
##              plan_self 1.27 [1.22,  1.32]         1.13      0.79
##           plan_partner 1.28 [1.24,  1.34]         1.13      0.78
##                    day 1.04 [1.02,  1.09]         1.02      0.96
##       weartime_self_cw 1.04 [1.02,  1.09]         1.02      0.96
##       weartime_self_cb 1.26 [1.22,  1.31]         1.12      0.79
##  Tolerance 95% CI
##      [0.86, 0.93]
##      [0.86, 0.92]
##      [0.89, 0.95]
##      [0.87, 0.93]
##      [0.85, 0.91]
##      [0.81, 0.87]
##      [0.24, 0.26]
##      [0.22, 0.24]
##      [0.76, 0.82]
##      [0.75, 0.81]
##      [0.92, 0.98]
##      [0.92, 0.98]
##      [0.76, 0.82]
## 
## Moderate Correlation
## 
##                   Term  VIF    VIF 95% CI Increased SE Tolerance
##     persuasion_self_cb 9.10 [8.58,  9.64]         3.02      0.11
##  persuasion_partner_cb 9.75 [9.20, 10.33]         3.12      0.10
##        pushing_self_cb 6.24 [5.90,  6.61]         2.50      0.16
##     pushing_partner_cb 6.00 [5.68,  6.35]         2.45      0.17
##  Tolerance 95% CI
##      [0.10, 0.12]
##      [0.10, 0.11]
##      [0.15, 0.17]
##      [0.16, 0.18]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy         88%
##     lognormal          9%
##       weibull          3%
## 
## Predicted Distribution of Response
## 
##                Distribution Probability
##                   lognormal         72%
##                     tweedie         16%
##  neg. binomial (zero-infl.)          9%
## 
## Divergences:
## 0 of 40000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 40000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## DHARMa:testOutliers with type = binomial may have inflated Type I error rates for integer-valued distributions. To get a more exact result, it is recommended to re-run testOutliers with type = 'bootstrap'. See ?testOutliers for details

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 29, observations = 3337, p-value < 2.2e-16
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.001198681 0.005251723
## sample estimates:
## outlier frequency (expected: 0.00296074318249925 ) 
##                                        0.008690441
if (do_priorsense) {
  gc()
  priorsense::powerscale_sensitivity(pa_obj_log, variable = priorsense_vars)
  priorsense::powerscale_plot_dens(pa_obj_log, variable = priorsense_vars)
  priorsense::powerscale_plot_ecdf(pa_obj_log, variable = priorsense_vars)
  priorsense::powerscale_plot_quantities(pa_obj_log, variable = priorsense_vars)
}
# rope range for lognormal model
rope_factor <- sd(log(pa_obj_log$data$pa_obj))
rope_range_log = c(-0.1 * rope_factor, 0.1 * rope_factor)

summarize_brms(
  pa_obj_log, 
  stats_to_report = stats_to_report,
  rope_range = rope_range_log,
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack)
## Sampling priors, please wait...
## Warning in summarize_brms(pa_obj_log, stats_to_report = stats_to_report, :
## Coefficients were exponentiated. Double check if this was intended.
exp(Est.) SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
Intercept 111.19*** 6.10 [99.90, 123.80] 1.000 [0.94, 1.07] 0.000 >100 Overwhelming Evidence 1.000 7302 14610
Within-Person Effects
Daily persuasion experienced 1.03 0.02 [ 1.00, 1.06] 0.954 [0.94, 1.07] 0.995 0.023 Very Strong Evidence for Null 1.000 31745 31489
Daily persuasion utilized (partner’s view) 1.02 0.02 [ 0.99, 1.05] 0.866 [0.94, 1.07] 0.998 0.010 Very Strong Evidence for Null 1.000 36811 32211
Daily pressure experienced 0.95 0.03 [ 0.88, 1.01] 0.950 [0.94, 1.07] 0.599 0.019 Very Strong Evidence for Null 1.000 58948 32857
Daily pressure utilized (partner’s view) 0.98 0.03 [ 0.92, 1.05] 0.692 [0.94, 1.07] 0.911 0.005 Very Strong Evidence for Null 1.000 62674 30987
Daily pushing experienced 1.01 0.02 [ 0.96, 1.07] 0.730 [0.94, 1.07] 0.976 0.007 Very Strong Evidence for Null 1.000 43490 31579
Daily pushing utilized (partner’s view) 1.00 0.02 [ 0.96, 1.04] 0.517 [0.94, 1.07] 0.997 0.005 Very Strong Evidence for Null 1.000 59076 31365
Day 0.97 0.03 [ 0.91, 1.04] 0.786 [0.94, 1.07] 0.842 0.006 Very Strong Evidence for Null 1.000 86932 31131
Own Actionplan 1.06* 0.03 [ 1.01, 1.12] 0.994 [0.94, 1.07] 0.527 0.103 Moderate Evidence for Null 1.000 84768 29763
Partner Actionplan 1.05 0.03 [ 1.00, 1.10] 0.972 [0.94, 1.07] 0.751 0.025 Very Strong Evidence for Null 1.000 80571 29234
Daily weartime 1.00*** 0.00 [ 1.00, 1.00] 1.000 [0.94, 1.07] 1.000 15.508 Strong Evidence 1.000 91257 29645
Between-Person Effects
Mean persuasion experienced 1.10 0.15 [ 0.84, 1.46] 0.761 [0.94, 1.07] 0.281 0.017 Very Strong Evidence for Null 1.001 8040 15329
Mean persuasion utilized (partner’s view) 0.98 0.14 [ 0.74, 1.30] 0.560 [0.94, 1.07] 0.352 0.013 Very Strong Evidence for Null 1.001 8085 14898
Mean pressure experienced 0.99 0.14 [ 0.74, 1.32] 0.530 [0.94, 1.07] 0.345 0.008 Very Strong Evidence for Null 1.001 10570 20107
Mean pressure utilized (partner’s view) 0.97 0.14 [ 0.74, 1.29] 0.577 [0.94, 1.07] 0.352 0.008 Very Strong Evidence for Null 1.000 9726 16655
Mean pushing experienced 0.95 0.19 [ 0.64, 1.44] 0.600 [0.94, 1.07] 0.246 0.011 Very Strong Evidence for Null 1.000 9570 17205
Mean pushing utilized (partner’s view) 1.22 0.24 [ 0.82, 1.83] 0.842 [0.94, 1.07] 0.154 0.018 Very Strong Evidence for Null 1.000 9223 16339
Mean weartime 1.00 0.00 [ 1.00, 1.00] 0.916 [0.94, 1.07] 1.000 0.017 Very Strong Evidence for Null 1.000 43760 32719
Random Effects
sd(Intercept) 0.29 0.04 [0.23, 0.39] NA NA NA NA NA 1.000 9588 17646
sd(Daily persuasion experienced) 0.05 0.01 [0.02, 0.08] NA NA NA NA NA 1.000 21159 18046
sd(Daily persuasion utilized (partner’s view)) 0.05 0.02 [0.02, 0.09] NA NA NA NA NA 1.000 17608 14641
sd(Daily pressure experienced) 0.04 0.03 [0.00, 0.14] NA NA NA NA NA 1.000 22275 23013
sd(Daily pressure utilized (partner’s view)) 0.03 0.03 [0.00, 0.11] NA NA NA NA NA 1.000 24379 20553
sd(Daily pushing experienced) 0.07 0.04 [0.01, 0.15] NA NA NA NA NA 1.000 10463 14026
sd(Daily pushing utilized (partner’s view)) 0.03 0.03 [0.00, 0.10] NA NA NA NA NA 1.000 15655 21319
Additional Parameters
sigma 0.57 0.01 [0.56, 0.59] NA NA NA NA NA 1.001 74585 28663 
plot(
  bayestestR::p_direction(pa_obj_log),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-3, 3)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(pa_obj_log, range = rope_range_log, ci = 1)
) + theme_bw()
## Possible multicollinearity between b_persuasion_partner_cb and
##   b_persuasion_self_cb (r = 0.89), b_pressure_self_cb and
##   b_persuasion_self_cb (r = 0.73), b_pressure_partner_cb and
##   b_persuasion_self_cb (r = 0.74), b_pressure_self_cb and
##   b_persuasion_partner_cb (r = 0.71), b_pressure_partner_cb and
##   b_persuasion_partner_cb (r = 0.77), b_pushing_partner_cb and
##   b_pushing_self_cb (r = 0.83). This might lead to inappropriate results.
##   See 'Details' in '?rope'.

# Nothing significant, no plots

Affect

range(df_double$aff, na.rm = T)
## [1] 0 5
hist(df_double$aff, breaks = 15)

Gaussian

formula <- bf(
  aff ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    plan_self + plan_partner +
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
)


prior1 <- c(
  brms::set_prior("normal(0, 5)", class = "b")
  ,brms::set_prior("normal(0, 20)", class = "Intercept", lb=1, ub=6)
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
  , brms::set_prior("student_t(3, 0, 2.5)", class = "sigma", lb = 0)
)

#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = gaussian()
#  )

#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

mood_gauss <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = gaussian(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("mood_gauss", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
mood_gauss_digest <- digest::digest(mood_gauss)
if (check_models) {
  check_brms(mood_gauss, log_pp_check = FALSE)
  DHARMa.check_brms.all(mood_gauss, integer = FALSE)
}
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                   Term  VIF     VIF 95% CI Increased SE Tolerance
##     persuasion_self_cw 1.12 [ 1.09,  1.17]         1.06      0.89
##  persuasion_partner_cw 1.08 [ 1.05,  1.12]         1.04      0.93
##       pressure_self_cw 1.09 [ 1.07,  1.14]         1.05      0.91
##    pressure_partner_cw 1.08 [ 1.05,  1.12]         1.04      0.93
##        pushing_self_cw 1.16 [ 1.12,  1.20]         1.08      0.86
##     pushing_partner_cw 1.15 [ 1.11,  1.19]         1.07      0.87
##              plan_self 1.28 [ 1.24,  1.33]         1.13      0.78
##           plan_partner 1.27 [ 1.23,  1.32]         1.13      0.79
##                    day 1.05 [ 1.02,  1.09]         1.02      0.96
##  Tolerance 95% CI
##      [0.86, 0.92]
##      [0.89, 0.95]
##      [0.88, 0.94]
##      [0.89, 0.95]
##      [0.83, 0.89]
##      [0.84, 0.90]
##      [0.75, 0.81]
##      [0.76, 0.81]
##      [0.92, 0.98]
## 
## Moderate Correlation
## 
##                 Term  VIF     VIF 95% CI Increased SE Tolerance
##     pressure_self_cb 5.89 [ 5.58,  6.22]         2.43      0.17
##  pressure_partner_cb 5.90 [ 5.59,  6.22]         2.43      0.17
##      pushing_self_cb 8.80 [ 8.33,  9.30]         2.97      0.11
##   pushing_partner_cb 8.82 [ 8.34,  9.32]         2.97      0.11
##  Tolerance 95% CI
##      [0.16, 0.18]
##      [0.16, 0.18]
##      [0.11, 0.12]
##      [0.11, 0.12]
## 
## High Correlation
## 
##                   Term   VIF     VIF 95% CI Increased SE Tolerance
##     persuasion_self_cb 14.99 [14.16, 15.87]         3.87      0.07
##  persuasion_partner_cb 14.88 [14.06, 15.75]         3.86      0.07
##  Tolerance 95% CI
##      [0.06, 0.07]
##      [0.06, 0.07]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##        cauchy         44%
##        normal         41%
##   exponential          6%
## 
## Predicted Distribution of Response
## 
##   Distribution Probability
##        tweedie         41%
##  beta-binomial         38%
##    half-cauchy          6%
## 
## Divergences:
## 0 of 40000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 40000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## 
##  DHARMa outlier test based on exact binomial test with approximate
##  expectations
## 
## data:  model.check
## outliers at both margin(s) = 29, observations = 3736, p-value =
## 1.68e-09
## alternative hypothesis: true probability of success is not equal to 0.001998002
## 95 percent confidence interval:
##  0.005204527 0.011129091
## sample estimates:
## frequency of outliers (expected: 0.001998001998002 ) 
##                                          0.007762313
if (do_priorsense) {
  gc()
  priorsense::powerscale_sensitivity(mood_gauss, variable = priorsense_vars)
  priorsense::powerscale_plot_dens(mood_gauss, variable = priorsense_vars)
  priorsense::powerscale_plot_ecdf(mood_gauss, variable = priorsense_vars)
  priorsense::powerscale_plot_quantities(mood_gauss, variable = priorsense_vars)
}
summarize_brms(
  mood_gauss, 
  stats_to_report = stats_to_report,
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = F) %>%
  print_df(rows_to_pack = rows_to_pack)
## Warning in compute_rope(range = rope_range): Collinearity detected. Some VIFs
## are > 10. This may invalidate ROPE inferences!

Check for Multicollinearity

Low Correlation

              Term  VIF     VIF 95% CI Increased SE Tolerance
persuasion_self_cw 1.12 [ 1.09,  1.17]         1.06      0.89

persuasion_partner_cw 1.08 [ 1.05, 1.12] 1.04 0.93 pressure_self_cw 1.09 [ 1.07, 1.14] 1.05 0.91 pressure_partner_cw 1.08 [ 1.05, 1.12] 1.04 0.93 pushing_self_cw 1.16 [ 1.12, 1.20] 1.08 0.86 pushing_partner_cw 1.15 [ 1.11, 1.19] 1.07 0.87 plan_self 1.28 [ 1.24, 1.33] 1.13 0.78 plan_partner 1.27 [ 1.23, 1.32] 1.13 0.79 day 1.05 [ 1.02, 1.09] 1.02 0.96 Tolerance 95% CI [0.86, 0.92] [0.89, 0.95] [0.88, 0.94] [0.89, 0.95] [0.83, 0.89] [0.84, 0.90] [0.75, 0.81] [0.76, 0.81] [0.92, 0.98]

Moderate Correlation

            Term  VIF     VIF 95% CI Increased SE Tolerance
pressure_self_cb 5.89 [ 5.58,  6.22]         2.43      0.17

pressure_partner_cb 5.90 [ 5.59, 6.22] 2.43 0.17 pushing_self_cb 8.80 [ 8.33, 9.30] 2.97 0.11 pushing_partner_cb 8.82 [ 8.34, 9.32] 2.97 0.11 Tolerance 95% CI [0.16, 0.18] [0.16, 0.18] [0.11, 0.12] [0.11, 0.12]

High Correlation

              Term   VIF     VIF 95% CI Increased SE Tolerance
persuasion_self_cb 14.99 [14.16, 15.87]         3.87      0.07

persuasion_partner_cb 14.88 [14.06, 15.75] 3.86 0.07 Tolerance 95% CI [0.06, 0.07] [0.06, 0.07]

## Sampling priors, please wait...
Est. SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
Intercept 3.67*** 0.11 [ 3.46, 3.88] 1.000 [-0.11, 0.11] 0.000 >100 Overwhelming Evidence 1.000 4390 8992
Within-Person Effects
Daily persuasion experienced 0.00 0.02 [-0.04, 0.04] 0.521 [-0.11, 0.11] 1.000 0.004 Very Strong Evidence for Null 1.000 39834 29740
Daily persuasion utilized (partner’s view) 0.02 0.02 [-0.02, 0.07] 0.813 [-0.11, 0.11] 1.000 0.006 Very Strong Evidence for Null 1.000 31489 29249
Daily pressure experienced -0.03 0.05 [-0.14, 0.07] 0.716 [-0.11, 0.11] 0.936 0.004 Very Strong Evidence for Null 1.000 44219 28749
Daily pressure utilized (partner’s view) -0.03 0.05 [-0.15, 0.08] 0.691 [-0.11, 0.11] 0.927 0.004 Very Strong Evidence for Null 1.000 37429 26049
Daily pushing experienced 0.01 0.03 [-0.06, 0.07] 0.575 [-0.11, 0.11] 0.999 0.004 Very Strong Evidence for Null 1.000 43555 30803
Daily pushing utilized (partner’s view) 0.07* 0.03 [ 0.00, 0.14] 0.977 [-0.11, 0.11] 0.896 0.033 Strong Evidence for Null 1.000 33119 28321
Day 0.26*** 0.06 [ 0.15, 0.37] 1.000 [-0.11, 0.11] 0.006 63.680 Very Strong Evidence 1.000 65116 29860
Own Actionplan 0.10** 0.04 [ 0.03, 0.18] 0.996 [-0.11, 0.11] 0.619 0.106 Moderate Evidence for Null 1.000 66702 29058
Partner Actionplan -0.03 0.04 [-0.11, 0.04] 0.791 [-0.11, 0.11] 0.983 0.004 Very Strong Evidence for Null 1.000 62577 30130
Daily weartime NA NA NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 0.34 0.28 [-0.23, 0.90] 0.887 [-0.11, 0.11] 0.150 0.029 Very Strong Evidence for Null 1.000 4803 9226
Mean persuasion utilized (partner’s view) 0.23 0.28 [-0.34, 0.79] 0.791 [-0.11, 0.11] 0.231 0.019 Very Strong Evidence for Null 1.000 4707 9230
Mean pressure experienced -0.30 0.27 [-0.85, 0.25] 0.863 [-0.11, 0.11] 0.182 0.016 Very Strong Evidence for Null 1.000 5644 11725
Mean pressure utilized (partner’s view) -0.32 0.27 [-0.87, 0.23] 0.876 [-0.11, 0.11] 0.172 0.017 Very Strong Evidence for Null 1.000 5485 11919
Mean pushing experienced 0.20 0.39 [-0.57, 0.99] 0.699 [-0.11, 0.11] 0.206 0.012 Very Strong Evidence for Null 1.002 6146 11770
Mean pushing utilized (partner’s view) 0.36 0.39 [-0.41, 1.15] 0.827 [-0.11, 0.11] 0.150 0.017 Very Strong Evidence for Null 1.002 6107 11814
Mean weartime NA NA NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.60 0.07 [0.47, 0.78] NA NA NA NA NA 1.001 7123 14047
sd(Daily persuasion experienced) 0.04 0.03 [0.00, 0.10] NA NA NA NA NA 1.000 10238 15606
sd(Daily persuasion utilized (partner’s view)) 0.07 0.03 [0.01, 0.13] NA NA NA NA NA 1.000 9815 7528
sd(Daily pressure experienced) 0.07 0.06 [0.00, 0.23] NA NA NA NA NA 1.000 16476 18018
sd(Daily pressure utilized (partner’s view)) 0.08 0.07 [0.00, 0.26] NA NA NA NA NA 1.000 14238 17562
sd(Daily pushing experienced) 0.05 0.04 [0.00, 0.14] NA NA NA NA NA 1.000 14637 15323
sd(Daily pushing utilized (partner’s view)) 0.07 0.05 [0.01, 0.16] NA NA NA NA NA 1.000 13455 13980
Additional Parameters
sigma 0.96 0.01 [0.94, 0.98] NA NA NA NA NA 1.000 61626 30661 
plot(
  bayestestR::p_direction(mood_gauss),
  priors = TRUE
)  + 
  coord_cartesian(xlim = c(-3, 3)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(mood_gauss, ci = 1)
) + theme_bw()
## Possible multicollinearity between b_pressure_self_cb and
##   b_persuasion_self_cb (r = 0.83), b_pressure_partner_cb and
##   b_persuasion_self_cb (r = 0.8), b_pressure_self_cb and
##   b_persuasion_partner_cb (r = 0.8), b_pressure_partner_cb and
##   b_persuasion_partner_cb (r = 0.83), b_pressure_partner_cb and
##   b_pressure_self_cb (r = 0.77), b_pushing_partner_cb and
##   b_pushing_self_cb (r = 0.89). This might lead to inappropriate results.
##   See 'Details' in '?rope'.

conditional_spaghetti(
  mood_gauss, 
  effects = c('pushing_partner_cw'),
  group_var = 'coupleID',
  plot_full_range = TRUE
)

$pushing_partner_cw

Reactance

range(df_double$reactance, na.rm = T) 
## [1] 0 5
hist(df_double$reactance, breaks = 7) 

hist(log(df_double$reactance+0.1), breaks = 10)

Ordinal

df_double$reactance_ordinal <- factor(df_double$reactance,
                                      levels = 0:5, 
                                      ordered = TRUE)

formula <- bf(
  reactance_ordinal ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    plan_self + plan_partner +
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
)


prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = cumulative() # HURDLE_CUMULATIVE
#  )


#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

reactance_ordinal <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::cumulative(),
  #control = list(adapt_delta = 0.95),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777
  , file = file.path("models_cache_brms", paste0("reactance_ordinal", suffix))
)
## Warning: Rows containing NAs were excluded from the model.
reactance_ordinal_digest <- digest::digest(reactance_ordinal)
if (check_models) {
  check_brms(reactance_ordinal)
  DHARMa.check_brms.all(reactance_ordinal, outliers_type = 'bootstrap')
}
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                   Term  VIF    VIF 95% CI Increased SE Tolerance
##       pressure_self_cw 4.55 [4.23,  4.90]         2.13      0.22
##    pressure_partner_cw 1.54 [1.46,  1.63]         1.24      0.65
##        pushing_self_cw 1.29 [1.23,  1.36]         1.13      0.78
##     pushing_partner_cw 1.09 [1.05,  1.15]         1.04      0.92
##     persuasion_self_cb 1.08 [1.05,  1.15]         1.04      0.92
##  persuasion_partner_cb 1.07 [1.04,  1.14]         1.03      0.93
##       pressure_self_cb 1.39 [1.32,  1.47]         1.18      0.72
##    pressure_partner_cb 1.19 [1.14,  1.26]         1.09      0.84
##                    day 1.98 [1.86,  2.11]         1.41      0.51
##  Tolerance 95% CI
##      [0.20, 0.24]
##      [0.61, 0.69]
##      [0.74, 0.81]
##      [0.87, 0.95]
##      [0.87, 0.96]
##      [0.88, 0.97]
##      [0.68, 0.76]
##      [0.79, 0.88]
##      [0.47, 0.54]
## 
## Moderate Correlation
## 
##                Term  VIF    VIF 95% CI Increased SE Tolerance Tolerance 95% CI
##  persuasion_self_cw 8.66 [8.02,  9.37]         2.94      0.12     [0.11, 0.12]
##     pushing_self_cb 5.86 [5.44,  6.33]         2.42      0.17     [0.16, 0.18]
##  pushing_partner_cb 7.32 [6.78,  7.91]         2.70      0.14     [0.13, 0.15]
##           plan_self 5.46 [5.06,  5.89]         2.34      0.18     [0.17, 0.20]
##        plan_partner 6.22 [5.76,  6.72]         2.49      0.16     [0.15, 0.17]
## 
## High Correlation
## 
##                   Term   VIF    VIF 95% CI Increased SE Tolerance
##  persuasion_partner_cw 10.18 [9.41, 11.02]         3.19      0.10
##  Tolerance 95% CI
##      [0.09, 0.11]
## Error in h(simpleError(msg, call)) : 
##   error in evaluating the argument 'x' in selecting a method for function 'print': Predictive errors are not defined for ordinal or categorical models.
## 
## Divergences:
## 3 of 40000 iterations ended with a divergence (0.0075%).
## Try increasing 'adapt_delta' to remove the divergences.
## 
## Tree depth:
## 0 of 40000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 6 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.

## 
##  DHARMa bootstrapped outlier test
## 
## data:  model.check
## outliers at both margin(s) = 0, observations = 756, p-value = 1
## alternative hypothesis: two.sided
##  percent confidence interval:
##  0.000000000 0.002645503
## sample estimates:
## outlier frequency (expected: 0.00044973544973545 ) 
##                                                  0
if (do_priorsense) {
  gc()
  priorsense::powerscale_sensitivity(reactance_ordinal, variable = priorsense_vars)
  priorsense::powerscale_plot_dens(reactance_ordinal, variable = priorsense_vars)
  priorsense::powerscale_plot_ecdf(reactance_ordinal, variable = priorsense_vars)
  priorsense::powerscale_plot_quantities(reactance_ordinal, variable = priorsense_vars)
}
summarize_brms(
  reactance_ordinal, 
  stats_to_report = stats_to_report,
  rope_range = c(-0.18, 0.18),
  model_rows_fixed = model_rows_fixed_ordinal,
  model_rows_random = model_rows_random_ordinal,
  model_rownames_fixed = model_rownames_fixed_ordinal,
  model_rownames_random = model_rownames_random_ordinal,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack_ordinal)
## Warning: There were 3 divergent transitions after warmup. Increasing
## adapt_delta above 0.8 may help. See
## http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
## Warning in compute_rope(range = rope_range): Collinearity detected. Some VIFs
## are > 10. This may invalidate ROPE inferences!

Check for Multicollinearity

Low Correlation

              Term  VIF    VIF 95% CI Increased SE Tolerance
  pressure_self_cw 4.55 [4.23,  4.90]         2.13      0.22

pressure_partner_cw 1.54 [1.46, 1.63] 1.24 0.65 pushing_self_cw 1.29 [1.23, 1.36] 1.13 0.78 pushing_partner_cw 1.09 [1.05, 1.15] 1.04 0.92 persuasion_self_cb 1.08 [1.05, 1.15] 1.04 0.92 persuasion_partner_cb 1.07 [1.04, 1.14] 1.03 0.93 pressure_self_cb 1.39 [1.32, 1.47] 1.18 0.72 pressure_partner_cb 1.19 [1.14, 1.26] 1.09 0.84 day 1.98 [1.86, 2.11] 1.41 0.51 Tolerance 95% CI [0.20, 0.24] [0.61, 0.69] [0.74, 0.81] [0.87, 0.95] [0.87, 0.96] [0.88, 0.97] [0.68, 0.76] [0.79, 0.88] [0.47, 0.54]

Moderate Correlation

           Term  VIF    VIF 95% CI Increased SE Tolerance Tolerance 95% CI

persuasion_self_cw 8.66 [8.02, 9.37] 2.94 0.12 [0.11, 0.12] pushing_self_cb 5.86 [5.44, 6.33] 2.42 0.17 [0.16, 0.18] pushing_partner_cb 7.32 [6.78, 7.91] 2.70 0.14 [0.13, 0.15] plan_self 5.46 [5.06, 5.89] 2.34 0.18 [0.17, 0.20] plan_partner 6.22 [5.76, 6.72] 2.49 0.16 [0.15, 0.17]

High Correlation

              Term   VIF    VIF 95% CI Increased SE Tolerance

persuasion_partner_cw 10.18 [9.41, 11.02] 3.19 0.10 Tolerance 95% CI [0.09, 0.11]

## Sampling priors, please wait...
OR SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
Intercept NA NA NA NA NA NA NA NA NA NA NA
Intercept[1] 3.41*** 1.06 [ 1.88, 6.37] 1.000 [0.84, 1.20] 0.000 >100 Overwhelming Evidence 1.000 29678 29857
Intercept[2] 7.41*** 2.37 [ 4.00, 14.12] 1.000 [0.84, 1.20] 0.000 >100 Overwhelming Evidence 1.000 29774 29014
Intercept[3] 20.74*** 7.07 [ 10.82, 41.38] 1.000 [0.84, 1.20] 0.000 >100 Overwhelming Evidence 1.000 30595 29058
Intercept[4] 90.95*** 35.54 [ 43.19, 202.38] 1.000 [0.84, 1.20] 0.000 >100 Overwhelming Evidence 1.000 32305 29085
Intercept[5] 3085.94*** 2038.94 [928.01, 13025.16] 1.000 [0.84, 1.20] 0.000 >100 Overwhelming Evidence 1.000 41858 29149
Within-Person Effects
Daily persuasion experienced 0.84* 0.07 [ 0.71, 0.99] 0.980 [0.84, 1.20] 0.555 0.390 Weak Evidence for Null 1.000 35488 27412
Daily persuasion utilized (partner’s view) 1.02 0.10 [ 0.83, 1.23] 0.585 [0.84, 1.20] 0.923 0.045 Strong Evidence for Null 1.000 32417 27395
Daily pressure experienced 1.84* 0.36 [ 1.16, 2.66] 0.992 [0.84, 1.20] 0.030 2.005 Weak Evidence 1.000 19813 21668
Daily pressure utilized (partner’s view) 1.24 0.30 [ 0.69, 2.11] 0.806 [0.84, 1.20] 0.374 0.086 Strong Evidence for Null 1.000 20861 17900
Daily pushing experienced 1.21 0.13 [ 0.98, 1.53] 0.962 [0.84, 1.20] 0.449 0.211 Moderate Evidence for Null 1.000 28732 28538
Daily pushing utilized (partner’s view) 0.94 0.12 [ 0.71, 1.23] 0.687 [0.84, 1.20] 0.773 0.049 Strong Evidence for Null 1.000 29004 22972
Day 1.49 0.52 [ 0.75, 2.95] 0.874 [0.84, 1.20] 0.214 0.074 Strong Evidence for Null 1.000 42893 30903
Own Actionplan 0.85 0.24 [ 0.48, 1.49] 0.718 [0.84, 1.20] 0.406 0.047 Strong Evidence for Null 1.000 42179 29531
Partner Actionplan 0.92 0.24 [ 0.56, 1.54] 0.619 [0.84, 1.20] 0.493 0.041 Strong Evidence for Null 1.000 43578 29710
Daily weartime NA NA NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 1.08 0.58 [ 0.36, 3.23] 0.560 [0.84, 1.20] 0.261 0.057 Strong Evidence for Null 1.000 19038 26032
Mean persuasion utilized (partner’s view) 1.47 0.92 [ 0.45, 5.24] 0.736 [0.84, 1.20] 0.189 0.068 Strong Evidence for Null 1.000 20617 27126
Mean pressure experienced 3.74* 2.15 [ 1.23, 12.07] 0.989 [0.84, 1.20] 0.017 0.783 Weak Evidence for Null 1.000 17768 25613
Mean pressure utilized (partner’s view) 1.14 0.69 [ 0.33, 3.70] 0.582 [0.84, 1.20] 0.224 0.053 Strong Evidence for Null 1.000 18385 24862
Mean pushing experienced 1.17 0.92 [ 0.25, 5.98] 0.580 [0.84, 1.20] 0.176 0.056 Strong Evidence for Null 1.000 15789 21244
Mean pushing utilized (partner’s view) 0.08* 0.08 [ 0.01, 0.55] 0.994 [0.84, 1.20] 0.006 1.752 Weak Evidence 1.000 18548 24773
Mean weartime NA NA NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 0.82 0.20 [0.48, 1.30] NA NA NA NA NA 1.000 14138 22006
sd(Daily persuasion experienced) 0.17 0.12 [0.01, 0.43] NA NA NA NA NA 1.000 7515 11752
sd(Daily persuasion utilized (partner’s view)) 0.21 0.14 [0.01, 0.52] NA NA NA NA NA 1.000 8586 11983
sd(Daily pressure experienced) 0.56 0.25 [0.10, 1.18] NA NA NA NA NA 1.000 8771 6371
sd(Daily pressure utilized (partner’s view)) 0.44 0.40 [0.02, 1.59] NA NA NA NA NA 1.000 8991 16845
sd(Daily pushing experienced) 0.21 0.14 [0.01, 0.51] NA NA NA NA NA 1.000 10556 12999
sd(Daily pushing utilized (partner’s view)) 0.15 0.14 [0.01, 0.63] NA NA NA NA NA 1.000 14325 18401
Additional Parameters
sigma NA NA NA NA NA NA NA NA NA NA NA
disc 1.00 0.00 [1.00, 1.00] NA NA NA NA NA NA NA NA 
plot(
  bayestestR::p_direction(reactance_ordinal),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-6, 6)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(reactance_ordinal, range = c(-0.18, 0.18), ci = 1)
) + theme_bw()
## Possible multicollinearity between b_Intercept[4] and b_Intercept[2] (r
##   = 0.81), b_Intercept[4] and b_Intercept[3] (r = 0.86),
##   b_pressure_self_cb and b_persuasion_self_cb (r = 0.73),
##   b_pressure_partner_cb and b_persuasion_partner_cb (r = 0.81). This might
##   lead to inappropriate results. See 'Details' in '?rope'.

conditional_spaghetti(
  reactance_ordinal, 
  effects = c('persuasion_self_cw', 'pressure_self_cw')
  , group_var = 'coupleID'
  #, n_groups = 15
  , plot_full_range = T
)

\(persuasion_self_cw <img src="01_FinalModelsPlan_files/figure-html/report_reactance_ordinal-3.png" width="2400" />\)pressure_self_cw

Binary

introduce_binary_reactance <- function(data) {
  data$is_reactance <- factor(data$reactance > 0, levels = c(FALSE, TRUE), labels = c(0, 1))
  return(data)
}



df_double <- introduce_binary_reactance(df_double)
if (use_mi) {
  for (i in seq_along(implist)) {
    implist[[i]] <- introduce_binary_reactance(implist[[i]])
  }
}


formula <- bf(
  is_reactance ~ 
    persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw +
    
    persuasion_self_cb + persuasion_partner_cb +
    pressure_self_cb + pressure_partner_cb +
    pushing_self_cb + pushing_partner_cb +
    
    plan_self + plan_partner +
    day +
    
    # Random effects
    (persuasion_self_cw + persuasion_partner_cw +
    pressure_self_cw + pressure_partner_cw +
    pushing_self_cw + pushing_partner_cw | coupleID)
  
  , decomp = 'QR'
  #, autocor = autocor_str
  )



prior1 <- c(
  brms::set_prior("normal(0, 2.5)", class = "b")
  , brms::set_prior("normal(0, 10)", class = "Intercept", lb=0, ub=5) 
  , brms::set_prior("normal(0, 2)", class = "sd", group = "coupleID", lb = 0)
)


#brms::validate_prior(
#  prior1, 
#  formula = formula, 
#  data = df_double, 
#  family = bernoulli()
#  )



#df_minimal <- df_double[, c("AorB", all.vars(as.formula(formula)))]

is_reactance <- my_brm(
  mi = use_mi, 
  imputed_data = implist,
  
  formula = formula, 
  prior = prior1,
  data = df_double, 
  family = brms::bernoulli(),
  #control = list(adapt_delta = 0.95, max_treedepth = 15),
  iter = iterations,
  warmup = warmup,
  chains = 4,
  cores = 4,
  seed = 7777,
  file = file.path("models_cache_brms", paste0("is_reactance", suffix))
  #, file_refit = 'always'
)
## Warning: Rows containing NAs were excluded from the model.
is_reactance_digest <- digest::digest(is_reactance)
if (check_models) {
  check_brms(is_reactance)
  DHARMa.check_brms.all(is_reactance, integer = FALSE)
}
## # Check for Multicollinearity
## 
## Low Correlation
## 
##                   Term  VIF   VIF 95% CI Increased SE Tolerance
##     persuasion_self_cw 1.06 [1.03, 1.13]         1.03      0.94
##  persuasion_partner_cw 1.15 [1.10, 1.22]         1.07      0.87
##       pressure_self_cw 1.05 [1.02, 1.12]         1.02      0.95
##    pressure_partner_cw 1.06 [1.02, 1.13]         1.03      0.95
##        pushing_self_cw 1.26 [1.20, 1.33]         1.12      0.79
##     pushing_partner_cw 1.15 [1.11, 1.22]         1.07      0.87
##     persuasion_self_cb 2.70 [2.52, 2.89]         1.64      0.37
##  persuasion_partner_cb 3.34 [3.12, 3.59]         1.83      0.30
##       pressure_self_cb 2.97 [2.78, 3.19]         1.72      0.34
##    pressure_partner_cb 3.14 [2.93, 3.37]         1.77      0.32
##        pushing_self_cb 2.20 [2.06, 2.35]         1.48      0.46
##     pushing_partner_cb 2.05 [1.93, 2.19]         1.43      0.49
##              plan_self 1.67 [1.58, 1.78]         1.29      0.60
##           plan_partner 1.55 [1.46, 1.64]         1.24      0.65
##                    day 1.06 [1.03, 1.13]         1.03      0.95
##  Tolerance 95% CI
##      [0.88, 0.97]
##      [0.82, 0.91]
##      [0.89, 0.98]
##      [0.89, 0.98]
##      [0.75, 0.83]
##      [0.82, 0.90]
##      [0.35, 0.40]
##      [0.28, 0.32]
##      [0.31, 0.36]
##      [0.30, 0.34]
##      [0.43, 0.48]
##      [0.46, 0.52]
##      [0.56, 0.63]
##      [0.61, 0.68]
##      [0.89, 0.98]
## # Distribution of Model Family
## 
## Predicted Distribution of Residuals
## 
##  Distribution Probability
##       tweedie         34%
##        normal         19%
##     bernoulli          9%
## 
## Predicted Distribution of Response
## 
##   Distribution Probability
##      bernoulli         81%
##  beta-binomial         12%
##       binomial          6%
## 
## Divergences:
## 0 of 40000 iterations ended with a divergence.
## 
## Tree depth:
## 0 of 40000 iterations saturated the maximum tree depth of 10.
## 
## Energy:
## E-BFMI indicated no pathological behavior.

## Using 10 posterior draws for ppc type 'ecdf_overlay' by default.

## Using 10 posterior draws for ppc type 'dens_overlay' by default.

## Warning: Found 32 observations with a pareto_k > 0.7 in model 'model'. We
## recommend to set 'moment_match = TRUE' in order to perform moment matching for
## problematic observations.

## 
##  DHARMa outlier test based on exact binomial test with approximate
##  expectations
## 
## data:  model.check
## outliers at both margin(s) = 0, observations = 756, p-value = 0.4141
## alternative hypothesis: true probability of success is not equal to 0.001998002
## 95 percent confidence interval:
##  0.000000000 0.004867585
## sample estimates:
## frequency of outliers (expected: 0.001998001998002 ) 
##                                                    0
if (do_priorsense) {
  gc()
  priorsense::powerscale_sensitivity(is_reactance, variable = priorsense_vars)
  priorsense::powerscale_plot_dens(is_reactance, variable = priorsense_vars)
  priorsense::powerscale_plot_ecdf(is_reactance, variable = priorsense_vars)
  priorsense::powerscale_plot_quantities(is_reactance, variable = priorsense_vars)
}
summarize_brms(
  is_reactance, 
  stats_to_report = stats_to_report,
  model_rows_fixed = model_rows_fixed,
  model_rows_random = model_rows_random,
  model_rownames_fixed = model_rownames_fixed,
  model_rownames_random = model_rownames_random,
  exponentiate = T) %>%
  print_df(rows_to_pack = rows_to_pack)
## Sampling priors, please wait...
OR SE 95% CI pd ROPE inside ROPE BF BF_Evidence Rhat Bulk_ESS Tail_ESS
Intercept 0.34** 0.11 [0.17, 0.65] 0.999 [0.83, 1.20] 0.004 9.815 Moderate Evidence 1.000 43983 32996
Within-Person Effects
Daily persuasion experienced 0.84 0.08 [0.68, 1.01] 0.968 [0.83, 1.20] 0.524 0.287 Moderate Evidence for Null 1.000 39116 35830
Daily persuasion utilized (partner’s view) 1.12 0.16 [0.84, 1.53] 0.779 [0.83, 1.20] 0.663 0.088 Strong Evidence for Null 1.000 30621 31810
Daily pressure experienced 1.98* 0.64 [1.02, 4.58] 0.978 [0.83, 1.20] 0.053 0.966 Weak Evidence for Null 1.000 27810 23603
Daily pressure utilized (partner’s view) 1.41 0.62 [0.54, 4.37] 0.795 [0.83, 1.20] 0.231 0.137 Moderate Evidence for Null 1.000 27990 25348
Daily pushing experienced 1.34* 0.18 [1.04, 1.76] 0.988 [0.83, 1.20] 0.206 0.656 Weak Evidence for Null 1.000 39055 33035
Daily pushing utilized (partner’s view) 0.92 0.18 [0.61, 1.37] 0.671 [0.83, 1.20] 0.601 0.073 Strong Evidence for Null 1.000 40541 30682
Day 1.68 0.67 [0.77, 3.69] 0.904 [0.83, 1.20] 0.161 0.101 Moderate Evidence for Null 1.000 64852 30859
Own Actionplan 0.87 0.28 [0.46, 1.62] 0.669 [0.83, 1.20] 0.392 0.048 Strong Evidence for Null 1.000 68276 30117
Partner Actionplan 0.88 0.26 [0.49, 1.58] 0.669 [0.83, 1.20] 0.425 0.048 Strong Evidence for Null 1.000 73521 28927
Daily weartime NA NA NA NA NA NA NA NA NA NA NA
Between-Person Effects
Mean persuasion experienced 1.89 1.26 [0.51, 7.72] 0.832 [0.83, 1.20] 0.135 0.110 Moderate Evidence for Null 1.000 29689 32175
Mean persuasion utilized (partner’s view) 1.99 1.51 [0.47, 9.47] 0.825 [0.83, 1.20] 0.129 0.100 Moderate Evidence for Null 1.000 29786 32205
Mean pressure experienced 26.35** 32.93 [2.64, 404.79] 0.998 [0.83, 1.20] 0.002 6.716 Moderate Evidence 1.000 23979 25932
Mean pressure utilized (partner’s view) 1.77 2.30 [0.11, 19.86] 0.669 [0.83, 1.20] 0.095 0.124 Moderate Evidence for Null 1.000 22964 28932
Mean pushing experienced 0.61 0.69 [0.06, 5.70] 0.672 [0.83, 1.20] 0.115 0.089 Strong Evidence for Null 1.000 23562 29427
Mean pushing utilized (partner’s view) 0.05* 0.06 [0.00, 0.55] 0.992 [0.83, 1.20] 0.007 1.543 Weak Evidence 1.000 23173 30242
Mean weartime NA NA NA NA NA NA NA NA NA NA NA
Random Effects
sd(Intercept) 1.18 0.26 [0.75, 1.79] NA NA NA NA NA 1.000 17774 26197
sd(Daily persuasion experienced) 0.21 0.14 [0.01, 0.51] NA NA NA NA NA 1.000 10723 17026
sd(Daily persuasion utilized (partner’s view)) 0.49 0.21 [0.12, 0.99] NA NA NA NA NA 1.001 13264 12658
sd(Daily pressure experienced) 1.10 0.56 [0.13, 2.47] NA NA NA NA NA 1.001 9043 8898
sd(Daily pressure utilized (partner’s view)) 0.86 0.70 [0.05, 2.83] NA NA NA NA NA 1.000 15190 21082
sd(Daily pushing experienced) 0.25 0.16 [0.02, 0.61] NA NA NA NA NA 1.000 14769 18126
sd(Daily pushing utilized (partner’s view)) 0.27 0.24 [0.01, 0.99] NA NA NA NA NA 1.000 16875 20309
Additional Parameters
sigma NA NA NA NA NA NA NA NA NA NA NA 
plot(
  bayestestR::p_direction(is_reactance),
  priors = TRUE
) + 
  coord_cartesian(xlim = c(-6, 6)) +
  theme_bw()
## Warning in `==.default`(dens$Parameter, parameter): longer object length is not
## a multiple of shorter object length
## Warning in is.na(e1) | is.na(e2): longer object length is not a multiple of
## shorter object length

plot(
  bayestestR::rope(is_reactance, ci = 1)
) + theme_bw()

conditional_spaghetti(
  is_reactance, 
  effects = c('pressure_self_cw', 'pushing_self_cw'),
  group_var = 'coupleID',
  plot_full_range = TRUE
)

\(pressure_self_cw <img src="01_FinalModelsPlan_files/figure-html/report_is_reactance-3.png" width="2400" />\)pushing_self_cw 

hypothesis(is_reactance, "exp(pressure_self_cw) > exp(pushing_self_cw)")
## Hypothesis Tests for class b:
##                 Hypothesis Estimate Est.Error CI.Lower CI.Upper Evid.Ratio
## 1 (exp(pressure_sel... > 0     0.83      1.01     -0.3     2.56       6.26
##   Post.Prob Star
## 1      0.86     
## ---
## 'CI': 90%-CI for one-sided and 95%-CI for two-sided hypotheses.
## '*': For one-sided hypotheses, the posterior probability exceeds 95%;
## for two-sided hypotheses, the value tested against lies outside the 95%-CI.
## Posterior probabilities of point hypotheses assume equal prior probabilities.

Report All Models

all_models <- report_side_by_side(
  pa_sub,
  pa_obj_log,
  mood_gauss,
  #reactance_ordinal,
  is_reactance,
  
  stats_to_report = c('CI'),
  
  model_rows_random = model_rows_random,
  model_rows_fixed = model_rows_fixed,
  model_rownames_random = model_rownames_random,
  model_rownames_fixed = model_rownames_fixed
)

[1] “pa_sub”

## Warning in summarize_brms(model, exponentiate = exponentiate, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this was
## intended.

[1] “pa_obj_log”

## Warning in summarize_brms(model, exponentiate = exponentiate, stats_to_report =
## stats_to_report, : Coefficients were exponentiated. Double check if this was
## intended.

[1] “mood_gauss” [1] “is_reactance”

# pretty printing

summary_all_models <- all_models %>%
  print_df(rows_to_pack = rows_to_pack) %>%
  add_header_above(
    c(" ", "Subjective MVPA Hurdle Lognormal" = 4,  
      "Device-Based MVPA Log (Gaussian)" = 2, 
      "Mood Gaussian" = 2,
      #"Reactance Ordinal" = 2,
      "Reactance Dichotome" = 2
  )
)


export_xlsx(
  summary_all_models, 
  rows_to_pack = rows_to_pack,
  file.path("Output", paste0("AllModels", suffix, ".xlsx")), 
  merge_option = 'header', 
  simplify_2nd_row = TRUE,
  line_above_rows = c(1,2),
  line_below_rows = c(-1)
)
## 
## Attaching package: 'rvest'
## The following object is masked from 'package:readr':
## 
##     guess_encoding
print(summary_all_models)
Subjective MVPA Hurdle Lognormal
Device-Based MVPA Log (Gaussian)
Mood Gaussian
Reactance Dichotome
exp(Est.)_hu pa_sub 95% CI_hu pa_sub exp(Est.)_nonzero pa_sub 95% CI_nonzero pa_sub exp(Est.) pa_obj_log 95% CI pa_obj_log Est. mood_gauss 95% CI mood_gauss OR is_reactance 95% CI is_reactance
Intercept 0.26*** [ 0.19, 0.35] 37.84*** [33.00, 43.33] 111.19*** [99.90, 123.80] 3.67*** [ 3.46, 3.88] 0.34** [0.17, 0.65]
Within-Person Effects
Daily persuasion experienced 1.58*** [ 1.39, 1.83] 1.03 [ 0.98, 1.09] 1.03 [ 1.00, 1.06] 0.00 [-0.04, 0.04] 0.84 [0.68, 1.01]
Daily persuasion utilized (partner’s view) 1.34*** [ 1.19, 1.52] 1.03 [ 0.99, 1.08] 1.02 [ 0.99, 1.05] 0.02 [-0.02, 0.07] 1.12 [0.84, 1.53]
Daily pressure experienced 0.96 [ 0.69, 1.30] 0.91 [ 0.82, 1.00] 0.95 [ 0.88, 1.01] -0.03 [-0.14, 0.07] 1.98* [1.02, 4.58]
Daily pressure utilized (partner’s view) 1.49* [ 1.05, 2.34] 0.95 [ 0.86, 1.03] 0.98 [ 0.92, 1.05] -0.03 [-0.15, 0.08] 1.41 [0.54, 4.37]
Daily pushing experienced 0.95 [ 0.71, 1.30] 0.99 [ 0.92, 1.05] 1.01 [ 0.96, 1.07] 0.01 [-0.06, 0.07] 1.34* [1.04, 1.76]
Daily pushing utilized (partner’s view) 1.29* [ 1.05, 1.62] 0.96 [ 0.90, 1.02] 1.00 [ 0.96, 1.04] 0.07* [ 0.00, 0.14] 0.92 [0.61, 1.37]
Day 0.86 [ 0.64, 1.16] 0.99 [ 0.88, 1.12] 0.97 [ 0.91, 1.04] 0.26*** [ 0.15, 0.37] 1.68 [0.77, 3.69]
Own Actionplan 9.45*** [ 7.76, 11.57] 1.32*** [ 1.21, 1.45] 1.06* [ 1.01, 1.12] 0.10** [ 0.03, 0.18] 0.87 [0.46, 1.62]
Partner Actionplan 1.17 [ 0.96, 1.42] 1.08 [ 0.99, 1.17] 1.05 [ 1.00, 1.10] -0.03 [-0.11, 0.04] 0.88 [0.49, 1.58]
Daily weartime NA NA NA NA 1.00*** [ 1.00, 1.00] NA NA NA NA
Between-Person Effects
Mean persuasion experienced 1.49 [ 0.72, 3.06] 0.99 [ 0.74, 1.34] 1.10 [ 0.84, 1.46] 0.34 [-0.23, 0.90] 1.89 [0.51, 7.72]
Mean persuasion utilized (partner’s view) 1.40 [ 0.69, 2.85] 0.95 [ 0.71, 1.28] 0.98 [ 0.74, 1.30] 0.23 [-0.34, 0.79] 1.99 [0.47, 9.47]
Mean pressure experienced 0.39* [ 0.17, 0.88] 1.21 [ 0.86, 1.69] 0.99 [ 0.74, 1.32] -0.30 [-0.85, 0.25] 26.35** [2.64, 404.79]
Mean pressure utilized (partner’s view) 0.49 [ 0.21, 1.14] 0.94 [ 0.66, 1.33] 0.97 [ 0.74, 1.29] -0.32 [-0.87, 0.23] 1.77 [0.11, 19.86]
Mean pushing experienced 1.10 [ 0.38, 3.17] 1.24 [ 0.80, 1.94] 0.95 [ 0.64, 1.44] 0.20 [-0.57, 0.99] 0.61 [0.06, 5.70]
Mean pushing utilized (partner’s view) 1.96 [ 0.69, 5.73] 1.34 [ 0.85, 2.11] 1.22 [ 0.82, 1.83] 0.36 [-0.41, 1.15] 0.05* [0.00, 0.55]
Mean weartime NA NA NA NA 1.00 [ 1.00, 1.00] NA NA NA NA
Random Effects
sd(Intercept) 0.74 [0.57, 0.99] 0.30 [0.23, 0.40] 0.29 [0.23, 0.39] 0.60 [0.47, 0.78] 1.18 [0.75, 1.79]
sd(Daily persuasion experienced) 0.21 [0.03, 0.39] 0.11 [0.07, 0.17] 0.05 [0.02, 0.08] 0.04 [0.00, 0.10] 0.21 [0.01, 0.51]
sd(Daily persuasion utilized (partner’s view)) 0.15 [0.01, 0.34] 0.08 [0.04, 0.13] 0.05 [0.02, 0.09] 0.07 [0.01, 0.13] 0.49 [0.12, 0.99]
sd(Daily pressure experienced) 0.20 [0.01, 0.74] 0.07 [0.00, 0.23] 0.04 [0.00, 0.14] 0.07 [0.00, 0.23] 1.10 [0.13, 2.47]
sd(Daily pressure utilized (partner’s view)) 0.24 [0.01, 0.94] 0.06 [0.00, 0.18] 0.03 [0.00, 0.11] 0.08 [0.00, 0.26] 0.86 [0.05, 2.83]
sd(Daily pushing experienced) 0.59 [0.31, 1.01] 0.09 [0.01, 0.17] 0.07 [0.01, 0.15] 0.05 [0.00, 0.14] 0.25 [0.02, 0.61]
sd(Daily pushing utilized (partner’s view)) 0.18 [0.01, 0.50] 0.07 [0.01, 0.15] 0.03 [0.00, 0.10] 0.07 [0.01, 0.16] 0.27 [0.01, 0.99]
Additional Parameters
sigma NA NA 0.68 [0.65, 0.70] 0.57 [0.56, 0.59] 0.96 [0.94, 0.98] NA NA 
report::report_system()

Analyses were conducted using the R Statistical language (version 4.4.2; R Core Team, 2024) on Windows 11 x64 (build 26100)

report::report_packages()
  • beepr (version 2.0; Bååth R, 2024)
  • R.methodsS3 (version 1.8.2; Bengtsson H, 2003)
  • R.oo (version 1.27.0; Bengtsson H, 2003)
  • R.utils (version 2.12.3; Bengtsson H, 2023)
  • brms (version 2.22.0; Bürkner P, 2017)
  • digest (version 0.6.37; Eddelbuettel D, 2024)
  • Rcpp (version 1.0.13.1; Eddelbuettel D et al., 2024)
  • bayesplot (version 1.11.1; Gabry J, Mahr T, 2024)
  • lubridate (version 1.9.3; Grolemund G, Wickham H, 2011)
  • DHARMa (version 0.4.7; Hartig F, 2024)
  • wbCorr (version 0.1.22; Küng P, 2023)
  • see (version 0.9.0; Lüdecke D et al., 2021)
  • tibble (version 3.2.1; Müller K, Wickham H, 2023)
  • R (version 4.4.2; R Core Team, 2024)
  • openxlsx (version 4.2.7.1; Schauberger P, Walker A, 2024)
  • ggplot2 (version 3.5.1; Wickham H, 2016)
  • forcats (version 1.0.0; Wickham H, 2023)
  • stringr (version 1.5.1; Wickham H, 2023)
  • rvest (version 1.0.4; Wickham H, 2024)
  • tidyverse (version 2.0.0; Wickham H et al., 2019)
  • readxl (version 1.4.3; Wickham H, Bryan J, 2023)
  • dplyr (version 1.1.4; Wickham H et al., 2023)
  • purrr (version 1.0.2; Wickham H, Henry L, 2023)
  • readr (version 2.1.5; Wickham H et al., 2024)
  • xml2 (version 1.3.6; Wickham H et al., 2023)
  • tidyr (version 1.3.1; Wickham H et al., 2024)
  • knitr (version 1.49; Xie Y, 2024)
  • kableExtra (version 1.4.0; Zhu H, 2024)
report::cite_packages()